Modèle linéaire multiple

2eme-FA-EMS - BUT SD - E. Anakok

0.1 Problématique biologique

WarningObjectif

• Peut-on prédire le poids des poissons par leurs largeurs, leurs longueurs et leurs épaisseurs ?
• Tester le caractère significatif de la liaison : Y’a t’il un lien significatif entre le poids des poissons et leurs largeurs, leurs longueurs et leurs épaisseurs ?
• Y’a t’il un lien entre le poids du poisson et sa largeur après avoir pris en compte sa longueur et son epaisseur ?

• Données On a pour 20 brèmes péchées dans le lac Laengelmavesi en Finland leurs poids (en gramme) et leurs longeurs (en cm), leurs largeurs (en cm) leurs épaisseur (en cm).

1 Écriture et remarque sur le modèle

1.1 Écriture du modèle

1.1.0.1 Notations

On a \(n = 20\) observations. On note, pour \(1 \leq i \leq 20\)

• \(x_{1,i}\) la mesure de la longueur du poisson \(i\).
• \(x_{2,i}\) la mesure de la largeure du poisson \(i\).
• \(x_{3,i}\) la mesure de l’epaisseur du poisson \(i\).
• \(y_i\) la mesure du poids du poisson \(i\).

NoteDéfinition : Modèle de régression linéaire multiple

On suppose que \(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que: \[Y_i = \alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta + E_i,\quad 1 \leq i \leq n\]

• \(\alpha_1\) (resp \(\alpha_2\) , \(\alpha_3\) ) est un paramètre inconnu, l’effet de la longueur (resp largeur, épaisseur) sur le poids;
• \(\beta\) est un paramètre inconnu;
• \(E_i\) une variable aléatoire, \(E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\)
• Les \(x_i\) sont linéairement indépendants

1.2 Écriture plus générale du modèle

1.2.0.1 Notations

On a \(n = 20\) observations. On note, pour \(1 \leq i \leq 20\)

• \(x_{1}, \dots, x_{p}\), \(p\) variables explicatives et \(y\) une variable réponse
• \(y_i\) la mesure de la variable réponse de l’individu \(i\).
• \(\forall j \in \{ 1,\dots, p \}\), \(x_{j,i}\) la valeur de la variable \(x_j\) de l’individu \(i\)

NoteDéfinition : Modèle de régression linéaire multiple

On suppose que \(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que: \[Y_i = \sum_{j=1}^p \alpha_j x_{j,i} + \beta + E_i,\quad 1 \leq i \leq n\]

• \(\alpha_j\) est un paramètre inconnu, l’effet de la variable \(x_j\) sur le y;
• \(\beta\) est un paramètre inconnu
• \(E_i\) une variable aléatoire, \(E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\)
• Les \(x_i\) sont linéairement indépendants (en pratique vérifier avec le VIF)

1.3 Deux écritures équivalentes

CautionRemarques

Les deux modèles suivant sont équivalent :

• \(\forall i \in \{1, \dots, n\} \;\;\)

\[ Y_i = \sum_{j=1}^p \alpha_j x_{j,i} + \beta + E_i,\quad \ E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\]

• \(\forall i \in \{1, \dots, n\} \;\;\)

\[ Y_i \overset{i.i.d.}\sim \mathcal{N}\left(\sum_{j=1}^p \alpha_j x_{j,i} + \beta, \sigma^2\right)\]

1.4 Dans l’exemple des poissons

CautionRemarques

Les deux modèles suivant sont équivalent :

• \(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i = \alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta + E_i,\quad \ E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\]

• \(\forall i \in \{1, \dots, n\} \;\;\) \[Y_i \overset{i.i.d.}\sim \mathcal{N}\left(\alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta, \sigma^2\right)\]

1.5 Remarques sur les paramètres

CautionRemarques

Dans le modèle \(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i = \sum_{j=1}^p \alpha_j x_{j,i} + \beta + E_i,\quad\ E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\]

• \(\alpha_j\) représente l’accroissement de \(Y_i\) correspondant à l’accroissement d’une unité sur \(x_j\) si les autres variables explicatives sont fixées. Dans l’exemple des poissons : \(\alpha_1\) represente l’accroissement du poids correspondant à l’accroissement d’1 cm sur la longeurs du poisson si la largeur et l’épaisseur son fixé.

• Ce modèle est de dimension \(p + 1\). (Avec \(p+1\) paramètres d’espérance à estimer \(p\) pour \(\alpha\) et \(1\) pour\(\beta\))

• \(\forall i \in \{1, \dots, n\}\) \(Y_i\) se décompose en \(\color{red}{Y_i} = {\color{blue}{\underbrace{\displaystyle\sum_{j=1}^p \alpha_j x_{j,i} + \beta}_{{déterministe}{}}}} + \color{red}{\overbrace{E_i}^{aléatoire}}\)

2 Écriture matricielle du modèle

2.1 Rappel sur les matrices

1. Qu’est ce qu’une matrice ?

2. La multiplication matricielle

3. La transposée d’une matrice

4. L’inverse d’une matrice

5. Petits topo sur les vecteurs

2.2 Rappels d’algèbre

NoteRappel : Matrice

\(A = (a_{i,j})_{1\leq i \leq p, 1\leq j \leq q}\) est une matrice à coefficients réels de format \((p \times q)\)

\[A = \begin{bmatrix} a_{11} & a_{12} & \ldots & a_{1q}\\ a_{21} & a_{22} & \ldots & a_{2q}\\ \vdots & \vdots & & \vdots \\ a_{p1} & a_{p2} & \ldots & a_{pq} \end{bmatrix}\]

NoteRappel : Transposée d’une matrice

Si \(A = (a_{i,j})_{1\leq i \leq p, 1\leq j \leq q}\) est une matrice \((p \times q)\), on note alors \(A^t\) sa transposée, la matrice \((q \times p)\) définie par

\[A^t = \begin{bmatrix} a_{11} & a_{21} & \ldots & a_{p1}\\ a_{12} & a_{22} & \ldots & a_{p2}\\ \vdots & \vdots & & \vdots \\ a_{1q} & a_{2q} & \ldots & a_{pq} \end{bmatrix}\]

NoteRappel : Matrice symétrique

\(A\) est une matrice symétrique si :

• \(A\) est une matrice carré (\(p=q\))

• \(A^t = A\)

NoteRappel : Transposée du produit

Soit \(A\) et \(B\) deux matrices. On a

\[ (A B)^t = B^t A^t\]

2.3 Rappels d’algèbre : Vecteur

NoteRappel : Vecteurs

• On note V un vecteur à \(p\) éléments par

\[V = \begin{bmatrix} v_{1} \\ v_{2}\\ \vdots \\ v_{p} \end{bmatrix}= \left[{v_1, v_2, \ldots, v_p}\right]^t \]

2.4 Rappels d’algèbre : Vecteur

NoteRappel : Produit matrice-vecteur

Si \(V\) est un vecteur à \(q\) éléments et \(A\) une matrice \((p \times q)\) alors, \(U=AV\) est un vecteur à \(p\) éléments

\[ \begin{bmatrix} u_{1} \\ u_{2}\\ \vdots \\ u_{p} \end{bmatrix} = AV = \begin{bmatrix} a_{11} & a_{12} & \ldots & a_{1q}\\ a_{21} & a_{22} & \ldots & a_{2q}\\ \vdots & \vdots & & \vdots \\ a_{p1} & a_{p2} & \ldots & a_{pq} \end{bmatrix}\begin{bmatrix} v_{1} \\ v_{2}\\ \vdots \\ v_{p} \end{bmatrix} = \begin{bmatrix} \sum_{j=1}^{q} a_{1j} v_j \\ \sum_{j=1}^{q} a_{2j} v_j\\ \vdots \\ \sum_{j=1}^{q} a_{pj} \end{bmatrix}\]

2.5 Rappels d’algèbre : Identité et inverse

NoteRappel : matrice identité

La matrice identité \(I_n\) est la matrice carré à \(n\) lignes et \(n\) colonnes, dont les coefficients diagonaux valent \(1\) et les autres valent \(0\).

\[I_n = \begin{bmatrix} 1 & 0 & \ldots & 0\\ 0 & 1 & \ldots & 0\\ \vdots & \vdots & & \vdots \\ 0 & 0 & \ldots & 1 \end{bmatrix}\]

NoteRappel : Inverse d’une matrice

Si \(A = (a_{i,j})_{1\leq i \leq p, 1\leq j \leq p}\) est une matrice carrée \((p \times p)\) inversible, alors on note son inverse \(A^{-1}\) la matrice \((p \times p)\) telle que

\[A\,A^{-1} = A^{-1}\,A = I_p \]

NoteRappel : Transposée de l’inverse

Si \(A\) est une matrice inversible, alors

\[(A^{-1})^t=(A^t)^{-1}\]

2.6 Vecteurs aléatoires

NoteRappel : Vecteurs aléatoires

Soient \(V_1,V_2,\dots,V_n\), \(n\) variables aléatoires réelles, alors

\[V = \begin{bmatrix} V_1 \\ V_2 \\ \vdots \\ V_n \end{bmatrix} \]

NoteRappel : Espérance d’un vecteur aléatoire

\(V\) a pour espérance

\[\mathbb{E}[V] = \begin{bmatrix} \mathbb{E}[V_1] \\ \mathbb{E}[V_2] \\ \vdots \\ \mathbb{E}[V_n] \end{bmatrix}\]

NoteRappel : Matrice de variance-covariance

\(V\) a pour matrice de variance-covariance (symétrique)

\[\mathbb{V}[V] = \mathbb{E}\left[\left(V- \mathbb{E}[V] \right)\left(V- \mathbb{E}[V] \right)^t \right] =\]

\[\begin{bmatrix} \mathbb{V}[V_1] & cov(V_1,V_2) & \dots& Cov(V_1,V_n) \\ cov(V_2,V_1) & \mathbb{V}[V_2] & \dots& Cov(V_2,V_n) \\ \vdots & \vdots & \ddots & \vdots \\ cov(V_n,V_1) & cov(V_n,V_2) & \dots& \mathbb{V}[V_n]\end{bmatrix}\]

2.7 Propriétés

TipThéorème

Soient \(\color{blue}{A}\) et \(\color{blue}{B}\) des matrices à coefficients non aléatoires et soit \(\color{red}{V}\) un vecteur aléatoire, alors

\[\begin{align} \mathbb{E}[\color{blue}{A}\color{red}{V}\color{blue}{B}] &= \color{blue}{A} \mathbb{E}[\color{red}{V}]\color{blue}{B}\\ \mathbb{V}[\color{blue}{A}\color{red}{V}] & = \color{blue}{A} \mathbb{V}[\color{red}{V}] \color{blue}{A^t} \end{align}\]

NoteRappel : Vecteurs gaussiens

Si \(V_1, V_2, \ldots, V_n\) sont gaussiennes et indépendantes, alors \(V = [V_1, \ldots, V_n]^t\) est un vecteur gaussien.

Si \(V_1, V_2, \ldots, V_n \overset{i.i.d}\sim \mathcal{N}(0,\sigma ^2)\), alors

\[V \sim \mathcal{N}\left(\vec{0}_n,\,\sigma^2 I_n\right)\] où \(I_n\) est la matrice identité d’ordre \(n\).

2.8 Écriture matricielle du modèle

TipThéorème

\[\left\{\begin{matrix} Y_1 = \beta\ +\ \alpha_1\,x_{1,1}\ +\ \alpha_2\,x_{2,1}\ +\ \ldots\ +\ \alpha_{p}\,x_{p, 1}\ +\ E_1\\ Y_2 = \beta\ +\ \alpha_1\,x_{1, 2}\ +\ \alpha_2\,x_{2,2}\ +\ \ldots\ +\ \alpha_{p}\,x_{p, 2}\ +\ E_2\\ \vdots \\ Y_n = \beta\ +\ \alpha_1\,x_{1,n}\ +\ \alpha_2\,x_{2, n}\ +\ \ldots\ +\ \alpha_{p}\,x_{p,n}\ +\ E_n \end{matrix} \right.\]

peut se réécrire sous la forme

\[ Y = X \theta + E \]

avec

\[Y =\begin{bmatrix} Y_1\\ Y_2\\ \vdots \\ Y_n\end{bmatrix},\quad X = \begin{bmatrix} 1 & x_{1,1} & \ldots & x_{p, 1}\\ 1 & x_{1, 2} & \ldots & x_{p, 2}\\ \vdots & \vdots & & \vdots \\ 1 & x_{1, n} & \ldots & x_{p, n} \end{bmatrix},\quad \theta = \begin{bmatrix} \beta\\ \alpha_1\\ \alpha_2\\ \vdots \\ \alpha_{p}\end{bmatrix},\quad E = \begin{bmatrix}E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}\]

2.9 Écriture matricielle du modèle

TipPropriété

• \(X\) est de taille \(n\times (p+1)\)

• \(\theta\) est de dimension \((p+1)\times 1\) : vecteur des paramètres d’espérance

• \(E\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[E \sim \mathcal{N}\left(\vec{0}_n,\,\sigma^2 I_n\right) \]

• \(Y\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[ Y \sim \mathcal{N}\left(X\theta,\,\sigma^2 I_n\right)\]

CautionRemarques

Si les vecteurs colonnes de \(X\) sont linéairement indépendants, \(X\) est de rang \(p+1\) et \(X^t X\) est inversible.

2.10 Exercice : écrire le modèle linéaire univarié sous forme matricielle

2.11 Correction

Le modèle de régression simple \[Y_i=\alpha x_i+\beta+E_i\] s’écrit matriciellement \(Y=X\theta+E\) avec \[Y = \begin{bmatrix} Y_1\\ Y_2\\ \vdots \\ Y_n\end{bmatrix},\quad X = \begin{bmatrix} 1 & x_{1} \\ 1 & x_{2} \\ \vdots & \\ 1 & x_{n} \\ \end{bmatrix},\quad \theta\ =\ \begin{bmatrix}\beta\\ \alpha\end{bmatrix},\quad E\ =\ \begin{bmatrix} E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}\]

2.12 Écriture du modèle linéaire simple

TipPropriété

• \(X\) est de taille \((n\times 2)\)

• \(\theta\) est de dimension \((2 \times 1)\) : vecteur des paramètres d’espérance

• \(E\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[E \sim \mathcal{N}\left(\vec{0}_n,\,\sigma^2 I_n\right) \]

• \(Y\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[ Y \sim \mathcal{N}\left(X\theta,\,\sigma^2 I_n\right)\]

CautionRemarques

Si les \(x_i\) ne sont pas tous égaux, les deux vecteurs colonnes de \(X\) sont linéairement indépendants et \(X\) est de rang \(p = 2\).

2.13 Exemple sur les poissons

2.13.0.1 Notations

On a \(n = 20\) observations. On note, pour \(1 \leq i \leq 20\)

• \(x_{1,i}\) la mesure de la longueur du poisson \(i\).
• \(x_{2,i}\) la mesure de la largeure du poisson \(i\).
• \(x_{3,i}\) la mesure de l’epaisseur du poisson \(i\).
• \(y_i\) la mesure du poids du poisson \(i\).

NoteDéfinition : Modèle de régression linéaire multiple

On suppose que \(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que: \[Y_i = \alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta + E_i,\quad 1 \leq i \leq n\]

• \(\alpha_1\) (resp \(\alpha_2\) , \(\alpha_3\) ) est un paramètre inconnu, l’effet de la longueur (resp largeur, épaisseur) sur le poids;
• \(\beta\) est un paramètre inconnu;
• \(E_i\) une variable aléatoire, \(E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\)
• Les \(x_i\) sont linéairement indépendants

2.14 Comment écrire le modèle :

\[\left\{\begin{matrix} Y_1 = \beta\ +\ \alpha_1\,x_{1,1}\ +\ \alpha_2\,x_{2,1}\ +\ \alpha_{3}\,x_{3, 1}\ +\ E_1\\ Y_2 = \beta\ +\ \alpha_1\,x_{1, 2}\ +\ \alpha_2\,x_{2,2}\ +\ \alpha_{3}\,x_{3, 2}\ +\ E_2\\ \vdots \\ Y_n = \beta\ +\ \alpha_1\,x_{1,n}\ +\ \alpha_2\,x_{2, n}\ +\ \alpha_{3}\,x_{3,n}\ +\ E_n \end{matrix} \right.\]

peut se réécrire sous la forme

\[ Y = X \theta + E \]

avec

\[Y =\begin{bmatrix} Y_1\\ Y_2\\ \vdots \\ Y_n\end{bmatrix},\quad X = \begin{bmatrix} 1 & x_{1,1} & x_{2,1} & x_{3, 1}\\ 1 & x_{1, 2} & x_{2,1} & x_{3, 2}\\ \vdots & \vdots & & \vdots \\ 1 & x_{1, n} & x_{2,1} & x_{3, n} \end{bmatrix},\quad \theta = \begin{bmatrix} \beta\\ \alpha_1\\ \alpha_2\\ \alpha_{3}\end{bmatrix},\quad E = \begin{bmatrix}E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}\]

2.15 Avec les données :

\[X = \]

Intercept	Length1	Height	Width
1	23.2	11.5200	4.0200
1	24.0	12.4800	4.3056
1	23.9	12.3778	4.6961
1	26.3	12.7300	4.4555
1	26.5	12.4440	5.1340
1	26.8	13.6024	4.9274
1	26.8	14.1795	5.2785
1	27.6	12.6700	4.6900
1	27.6	14.0049	4.8438
1	28.5	14.2266	4.9594

3 Estimation des paramètres d’espérance du modèle

3.1 Estimateurs des moindres carrés

WarningObjectif

On note \(\theta= [\beta,\alpha_1,...,\alpha_p]^t\).

Comme dans la regression linéaire simple on cherche à minimiser :

\[J(\theta) = \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2\] Écriture matricielle ?

3.2 Rappels d’algèbre : norme

NoteRappel : norme

Soit \[V\ =\ \begin{bmatrix} v_{1} \\ v_{2}\\ \vdots \\ v_{n}\end{bmatrix} \ =\ [v_1, \ldots, v_n]^t\] un vecteur de \(\mathbb{R}^n\). On appelle norme de \(V\) le réel

\[\displaystyle V^t V \ =\ \sum_{i=1}^{n} v_i^2\ =\ ||V||^2 \]

CautionRemarque

\[V^t V \neq V V^t\]

\[V V^t = \begin{bmatrix} v_1^2 & v_1\,v_2 & \ldots & v_1\,v_n\\ v_2\,v_1 & v_2^2 & \ldots & v_2\,v_n\\ \vdots & \vdots & & \vdots \\ v_n\,v_1 & v_n\,v_2 & \ldots & v_n^2 \end{bmatrix}\]

3.3 Estimateurs des moindres carrés

WarningObjectif

On souhaite trouver \(\theta= (\beta,\alpha_1,...\alpha_p)^t\) qui minimise

\[J(\theta)=\sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 = (Y-X\theta)^t(Y-X\theta) = ||Y-X\theta||^2\]

TipThéorème

Le \(\theta\) optimal est la solution du système de \(p+1\) équations à \(p+1\) inconnues

\[ \left\{\begin{array}{ccc} \frac{\partial \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 }{\partial \beta} & = & 0 \\ % \frac{\partial \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 }{\partial \alpha_1} & = & 0\\ \vdots & & \\ \frac{\partial \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 }{\partial \alpha_p} & = & 0 \end{array}\right. \]

3.4 À la recherche de \(\widehat{\theta}\)

\[\widehat{\theta} = \underset{\theta\in\mathbb{R}^{p+1}}{argmin} (Y-X\theta)^t(Y-X\theta) =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\]

\(\widehat{\theta} = [B,A_1, \dots, A_p]^t\) où \(B\) est l’estimateur de \(\beta\) et \(A_j\) l’estimateur de \(\alpha_j\)

\[\begin{align} (Y-X\theta)^t(Y-X\theta) &= (Y^t-\theta^tX^t)(Y-X\theta)\\ &= Y^tY - Y^tX\theta -\theta^tX^t Y -\theta^tX^tX\theta \\ &= Y^tY - 2 Y^tX\theta -\theta^tX^tX\theta \end{align}\]

En dérivant par rapport à \(\theta\) on obtient : \[-2Y^tX + 2\theta^tX^tX \] Le système \(-2Y^tX + 2\theta^tX^tX=0\) est bien le même que les \(p+1\) equations de la slide précédente

3.5 À la recherche de \(\widehat{\theta}\)

\[\begin{align} -2Y^tX + 2\widehat{\theta}^tX^tX = 0 &\iff \widehat{\theta}^tX^tX = Y^tX \\ \iff X^tX \widehat{\theta} = X^tY & \iff \widehat{\theta} = (X^tX)^{-1} X^tY \end{align}\]

CautionRemarques

\(X\) est de rang \(p+1\) donc \(X^tX\) est inversible.

TipThéorème

\[\widehat{\theta} =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\] a pour solution

\[\widehat{\theta} = (X^tX)^{-1} X^tY\]

3.6 Estimateurs des paramètres d’espérance

WarningObjectif : Estimateur des moindres carrés de \(\theta\)

\[\widehat{\theta} =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\]

TipThéorème

• Si \(\hat{\theta}= [B, A_1, \ldots, A_{p}]^t\) est solution du système et que \((X^t\,X)\) est inversible, alors

\[\hat{\theta}=(X^t\,X)^{-1}\,X^t\,Y\]

Paramètres	Estimateurs	Estimations
\([\beta,\alpha_1,...\alpha_p]\)	\([B, A_1, \ldots, A_{p}]\)	\([b,a_1,\cdots,a_p]\)

4 Validation des hypothèses et estimation avec R

4.1 Les 4 graphiques

fish <- read.table(file = "Fish.csv",
                   sep =",", header = TRUE) %>% filter(Species == "Bream") 

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
par(mfrow=c(2,2))
plot(mod)

4.2 Multicolinéarité des \(x_i\)

On vérifie la corrélation entre les \(x_i\)

library(corrplot)

corrplot 0.95 loaded

M = cor(fish[c("Length1","Height","Width","Weight")])
corrplot.mixed(M)

4.3 Vérification avec le VIF

NoteDéfinition : VIF

\[ VIF(x_i) = \frac{1}{1-R^2(x_i)} \] où \(R^2(x_i)\) est le \(R^2\) du modèle linéaire dans lequel on explique \(x_i\) par toutes les autres covariables (\(x_j\))

4.4 Exemple du VIF sur les poissons :

Calcul ‘à la main’

mod_length<- lm(Length1 ~ Height + Width , data =  fish)
vif_length <- 1/ ( 1-summary(mod_length)$r.squared)

Avec la fonction VIF :

library(car)

Le chargement a nécessité le package : carData

Attachement du package : 'car'

L'objet suivant est masqué depuis 'package:dplyr':

    recode

L'objet suivant est masqué depuis 'package:purrr':

    some

vif(mod)

  Length1    Height     Width 
 8.952978 12.123683  7.451746 

En pratique : on retire les covariales qui ont un VIF supérieur à 10 (ou 5 si on veut être plus strict)

4.5 Estimations des paramètres sur notre exemple :

fish <- read.table(file = "Fish.csv",
                   sep =",", header = TRUE) %>% filter(Species == "Bream") %>% 
  slice(1:20)

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
summary(mod)

Call:
lm(formula = Weight ~ Length1 + Height + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-183.872  -17.044   -0.495   24.591  101.032 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1048.57     187.86  -5.582 4.13e-05 ***
Length1        18.96      13.26   1.430    0.172    
Height         48.50      28.31   1.713    0.106    
Width          66.72      51.73   1.290    0.215    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 63.03 on 16 degrees of freedom
Multiple R-squared:  0.808, Adjusted R-squared:  0.772 
F-statistic: 22.44 on 3 and 16 DF,  p-value: 5.624e-06

4.6 Sur notre exemple (sans la fonction lm)

1. On récupère les matrices X et Y :

X1 <- fish[, c("Length1",  "Height",  "Width")] %>% 
  as.matrix() 
X <- cbind(Intercept = 1,X1)
Y <- fish[,"Weight"] %>% 
  as.matrix()

4.7 Sur notre exemple (sans la fonction lm)

2. On cherche maintenant \(\widehat{\theta}=(X^t\,X)^{-1}\,X^t\,Y\).

Le produit matriciel \(AB\) s’obtient en faisant A %*% B, La transposée d’une matrice A s’obtient avec la fonction t(A) et son inverse avec la fonction solve(A).

• \(X^t\,X\) s’obtient avec :

t(X) %*% X

• \((X^t\,X)^{-1}\) s’obtient avec :

solve(t(X) %*% X)

On a finalement \(\widehat{\theta}\) :

solve(t(X) %*% X)%*% t(X)%*% Y

4.8 On retrouve bien les mêmes valeurs

coefficients(mod)

(Intercept)     Length1      Height       Width 
-1048.57013    18.95867    48.49646    66.71559 

solve(t(X) %*% X)%*% t(X)%*% Y

                 [,1]
Intercept -1048.57013
Length1      18.95867
Height       48.49646
Width        66.71559

5 Estimateur de la variance

5.1 Dans le modèle linéaire simple

NoteRappel

Dans le modèle \(M_2\) : \(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i =\beta + \alpha x_i + E_i,\quad E_i \overset{i.i.d.}{\sim}\mathcal{N}(0, \sigma^2)\]

L’estimateur de la variance résiduelle était :

\[ S^2 = \frac{1}{n-2}\sum_{i=1}^n(Y_i-\widehat{Y_i})^2=\frac{1}{n-2}\sum_{i=1}^n(Y_i-Ax_i-B)^2= \frac{SCR(M_2)}{n-2}, \]

où \(\widehat{Y_i}=Ax_i+B\), la prévision (aléatoire) par le modèle de régression linéaire associée à \(x_i\).

5.2 Le paramètre de variance résiduelle

WarningObjectif

Dans le modèle \(M_{p+1}\):

\(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i =\beta + \sum_{j=1}^p\alpha_j x_{j,i} + E_i, \quad E_i \overset{i.i.d.}{\sim}\mathcal{N}(0, \sigma^2)\]

il y’a un paramètre de variance \(\sigma ^2\) : la variance residuelle.

Quel est son estimateur ?

5.3 Estimation de la variance residuelle du modèle

NoteDéfinition : Prévision aléatoire par le modèle \(M_{p+1}\)

\[\widehat{Y_i}=B+A_1x_{i,1}+\cdots+A_px_{i,p}=B+\sum_{j=1}^p A_jx_{i,j}\]

NoteDéfinition : Résidus aléatoires du modèle \(M_{p+1}\)

\[E_i=Y_i-\widehat{Y_i}\]

NoteDéfinition : Somme des carrés résiduels

\[SCR(M_{p+1})=\sum_{i=1}^n E_i^2=\sum_{i=1}^n\left(Y_i-(B+A_1x_{i,1}+\cdots+A_px_{i,p})\right)^2\]

5.4 D’un point de vue matriciel

NoteDéfinition : Vecteur (aléatoire) des prévisions

\[\widehat{Y}\ =\ \begin{bmatrix} \widehat{Y_1}\\ \widehat{Y_2}\\ \vdots \\ \widehat{Y_n}\end{bmatrix}=\ \begin{bmatrix} B+\sum_{j=1}^p A_jx_{1,j}\\ B+\sum_{j=1}^p A_jx_{2,j}\\ \vdots \\ B+\sum_{j=1}^p A_jx_{n,j}\end{bmatrix}=X\hat{\theta}\]

NoteDéfinition : Vecteur (aléatoire) des résidus

\[E\ =\ \begin{bmatrix} E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}=\ \begin{bmatrix} Y_1-\widehat{Y_1}\\ Y_2-\widehat{Y_2}\\ \vdots \\ Y_n-\widehat{Y_n}\end{bmatrix}=Y-X\hat{\theta}\]

NoteDéfinition : Somme des carrés résiduels

\[SCR(M_{p+1})=\sum_{i=1}^n E_i^2= ||E||^2=||Y-X\hat{\theta}||^2\]

5.5 Estimateur de la variance résiduelle \(\sigma^2\)

NoteDéfinition : Estimateur de \(\sigma^2\)

\[\begin{align} {S^2}_{(M_{p+1})}&=\frac{SCR(M_{p+1})}{n-(p+1)}=\frac{1}{n-p-1}\sum_{i=1}^n E_i^2\\ &=\frac{1}{n-p-1}\sum_{i=1}^n\left(Y_i-(B+\sum_{j=1}^p A_jx_{i,j})\right)^2=\frac{||Y - X\,\hat{\theta}||^2}{n-p-1} \end{align}\]

NoteDéfinition : Estimation de \(\sigma^2\)

Réalisation \(s^2\) de \({S^2}_{(M_{p+1})}\) sur les données

\[s^2=\frac{1}{n-p-1}\sum_{i=1}^n e_i(M_{p+1})^2=\sum_{i=1}^n(y_i-(b+a_1x_{i,1}+\cdots+a_px_{i,p}))^2\] où \(e_i(M_{p+1})\) sont les résidus observés dans le modèle \(M_{p+1}\)

6 Propriétés et lois des estimateurs

6.1 Propriétés et lois de \(S^2\)

TipThéorème

Sous les hypothèses du modèle linéaire gaussien, \({S^2}_{(M_{p+1})}\) est un estimateur sans biais de \(\sigma^2\) et on a \[\frac{(n-p-1){S^2}_{(M_{p+1})}}{\sigma^2}=\frac{\sum_{i=1}^n(Y_i-B-\sum_{j=1}^pA_jx_{i,j})^2}{\sigma^2}\sim\chi^2(n-p-1)\]

De plus \({S^2}_{(M_{p+1})}\) est indépendant de \(\hat{\theta}\).

6.2 Propriétés et loi de \(\hat{\theta}\)

TipThéorème

\(\hat{\theta}\) est un vecteur gaussien d’espérance \(\theta\) et de matrice de variance-covariance \[Var(\hat{\theta})=\sigma^2(X^t\,X)^{-1} \quad \mbox{estimée par} \quad S^2(\hat{\theta})={S^2}_{(M_{p+1})}(X^t\,X)^{-1}\] On montre que \[\widehat{\theta}\sim{\cal N}(\theta\,,\,\sigma^2(X^t\,X)^{-1})\]

En particulier, pour tout \(k=1,\cdots,p+1\), si \(c_{kk}={(X^t\,X)^{-1}}_{kk}\) désigne le \(k\)-ème élément diagonal de \((X^t\,X)^{-1}\), et \(\widehat{\theta}_k\) désigne le \(k\)-ème élément du vecteur \(\hat{\theta}\), on a \[\widehat{\theta}_k\sim{\cal N}(\theta_k\,,\,\sigma^2c_{kk})\]

En remplaçant \(\sigma^2\) par son estimateur on a, \[\frac{(\widehat{\theta}_k-\theta_k)}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim \mathcal{T}{(n-p-1)}\]

6.3 Démonstrations

On rappelle que \[ \hat{\theta}=(X^t\,X)^{-1} X^t Y \]

\[ E(\hat{\theta})=E((X^t\,X)^{-1} X^t Y)=(X^t\,X)^{-1} X^t E(Y)=(X^t\,X)^{-1} X^t X\theta=\theta \]

\[\begin{align} {\rm{Var}}(\hat{\theta})&={\rm{Var}}((X^t\,X)^{-1} X^t Y)\\ &= (X^t\,X)^{-1} X^t {\rm{Var}}(Y) (X^t\,X)^{-1} X^t)^t\\ &= (X^t\,X)^{-1} X^t \sigma^2I_n ((X^t\,X)^{-1} X^t)^t\\ &= \sigma^2(X^t\,X)^{-1} X^t (X^t\,X)^{-1} X^t)^t\\ &= \sigma^2(X^t\,X)^{-1} X^t ((X^t)^t) ((X^t\,X)^{-1})^t \\ &= \sigma^2(X^t\,X)^{-1} X^t X ((X^t\,X)^{t})^{-1} \\ &= \sigma^2 (X^t\,(X^{t})^t)^{-1} \\ &= \sigma^2 (X^t\,X)^{-1} \\ \end{align}\]

6.4 Test et intervalle de confiance

À partir du résultat pour tout \(k=1,\cdots,p+1\)

\[\frac{(\hat{\theta}_k-\theta_k)}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim \mathcal{T}{(n-p-1)}\]

• On peut construire un intervalle de confiance de chaque paramètre de régression \(\theta_k\), où \(\theta_k\) désigne l’un des paramètres d’espérance (\(\beta\), ou \(\alpha_1\), …, ou \(\alpha_p\)).

• On peut construire le test de Student de la nullité de chaque coefficient de régression \(H_0:\theta_k=0\) contre \(H_1:\theta_k\neq 0\) à partir de la statistique de test \[T_n=\frac{\hat{\theta}_k}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim_{H_0} \mathcal{T}{(n-p-1)}\]

6.5 Intervalle de confiance de \(\theta_k\)

TipThéorème : Estimateurs par intervalle de niveau de confiance \(1-\delta\) des \(\theta_k\)

À partir des lois de \(\theta_k\) et \(S^2_{(M_{p+1})}\), on obtient:

\[IC_{1-\delta}(\theta_k) = \begin{align} \left[\hat{\theta}_k-t_{1-\frac{\delta}{2}} \sqrt{{S^2}_{(M_{p+1})}c_{kk}};\quad\hat{\theta}_k+t_{1-\frac{\delta}{2}} \sqrt{{S^2}_{(M_{p+1})}c_{kk}}\;\right]\\ \end{align}\]

où \(t_{1-\frac{\delta}{2}}\) est tel que \(\mathbb{P}\left(\mid \mathcal{T}(n-p-1)\mid \leq t_{1-\frac{\delta}{2}}\right)=1-\delta\).

\(t_{1-\frac{\delta}{2}}\) est le quantile d’ordre \(1-\frac{\delta}{2}\) de la loi \(\mathcal{T}(n-p-1)\).

6.6 Intervalle de confiance de \(\theta_k\) avec R :

confint(mod, level = 0.95)

                   2.5 %     97.5 %
(Intercept) -1446.805131 -650.33513
Length1        -9.153972   47.07131
Height        -11.528580  108.52149
Width         -42.941396  176.37258

7 V. Test dans le modèle de régression linéaire multiple

7.1 Test dans le modèle de régression multiple

WarningObjectifs

• Tester si un coefficient \(\theta_k\) est significativement différent de \(0\).
• Tester si la contribution globale des variables explicatives est significative pour expliquer \(Y\).
• Déterminer la contribution effective des variables explicatives dans la modélisation de \(Y\).
• Tester si un ensemble de \(q\) variables explicatives (\(q<p\)) ne suffit pas à modéliser correctement \(Y\).

7.2 Test de Student pour \(\theta_k\)

1. Statistique de test et loi sous \(H_0\)

\[{T_n= \frac{\hat{\theta}_k}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim_{H_0}\mathcal{T}(n- p -1)}\]

2. Zone de rejet

\[R_\alpha = \{T> t_{1-\frac{\delta}{2}} \} \]

On rejette \(H_0\) si \(t_{obs}> t_{1-\frac{\delta}{2}}\)

3. Application numérique

On calcul \(t_{obs} = \frac{a}{\sqrt{{s^2}_{(M_{p+1})}c_{kk}}}\) la réalisation de \(T_n\).

On compare avec \(t_{1-\frac{\delta}{2}}\) et on conclut.

7.3 \(p\)-valeur du test de Student.

TipCalcul de la \(p\)-valeur

\[p_c=\mathbb{P}_{H_0}(\mid T_n\mid >\mid t_n\mid)=2(1-\mathbb{P}(T_n\leq |t_n|))\] où \(T_n\sim \mathcal{T}(n-p -1)\)

NoteInteprétation

Pour un risque de 1ère espèce \(\delta\) fixé acceptable (par ex \(\delta=5\%\))

• si \(p_c <\delta\), le test de niveau \(\delta\) est significatif (liaison significative), on rejette \(H_0\).
• si \(p_c >\delta\), le test de niveau \(\delta\) n’est pas significatif (liaison non significative), on ne rejette pas \(H_0\).

7.4 Test de Student avec R

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
summary(mod)

Call:
lm(formula = Weight ~ Length1 + Height + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-183.872  -17.044   -0.495   24.591  101.032 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1048.57     187.86  -5.582 4.13e-05 ***
Length1        18.96      13.26   1.430    0.172    
Height         48.50      28.31   1.713    0.106    
Width          66.72      51.73   1.290    0.215    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 63.03 on 16 degrees of freedom
Multiple R-squared:  0.808, Adjusted R-squared:  0.772 
F-statistic: 22.44 on 3 and 16 DF,  p-value: 5.624e-06

7.5 Test de la contribution globale des variables explicatives

NoteGénéralisation du test de Fisher à \(p\) variables explicatives

\(H_0\) : aucune des \(p\) variables explicatives n’a d’influence sur \(Y\).

\(H_1\) : au moins une des variables explicatives contribue à expliquer \(Y\).

Autre formulation :

Test du modèle constant \[H_0\ :\ \text{modèle}\ M_1:Y_i=\beta+E_i,\quad E_i\ \overset{i.i.d.}{\sim}\ {\cal N}(0,\sigma^2)\] contre le modèle complet

\[H_1\ :\ \text{modèle}\ M_{p+1}:Y_i=\beta+\sum_{j=1}^p\alpha_j x_{i,j}+E_i,\quad E_i\ \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

Autre formulation :

Test de \(H_0:\ \forall j=1,\cdots,p,\ \alpha_j=0\quad\) contre \(\quad H_1:\exists j\ \alpha_j\neq 0\)

7.6 Test de Fisher de la contribution globale des variables explicatives

NoteDéfinition : \(SCT\)

La variabilité de \(Y\) sans tenir compte du modèle.

\[\color{purple}{SCT =\displaystyle\sum_{i = 1}^n( Y_i - \bar{Y})^2}\]

NoteDéfinition : \(SCM(M_{p+1})\)

Partie de la variabilité de \(Y\) expliquée par le modèle \(M_{p+1}\).

\[\color{blue}{SCM(M_{p+1}) = \displaystyle\sum_{i=1}^n(\widehat{Y_i}(M_{p+1})-\bar{Y})^2}\]

NoteDéfinition : \(SCR(M_{p+1})\)

Partie de la variabilité de \(Y\) qui n’est pas expliquée par le modèle \(M_{p+1}\).

\[\color{red}{\begin{align}SCR(M_{p+1}) &= \displaystyle\sum_{i=1}^n(Y_i-\widehat{Y_i}(M_{p+1}))^2\\&=\displaystyle\sum_{i=1}^n E_i(M_{p+1}) ^2\end{align}}\]

TipThéorème de décomposition de la variance du modèle \(M_{p+1}\)

\[\color{purple}{SCT}=\color{blue}{SCM(M_{p+1})}+ \color{red}{SCR(M_{p+1})}\]

7.7 Test de Fisher de la contribution globale des variables explicatives (suite)

CautionRemarques

Dans le modèle \(M_1\), \(\beta\) est estimée par \(\bar{Y}\) et \(\widehat{Y_i}(M_{1})=\bar{Y}\) donc \[SCR(M_1)=\sum_{i=1}^n(Y_i-\widehat{Y_i}(M_{1}))^2=\sum_{i=1}^n(Y_i-\bar{Y})^2=SCT\]

On a donc \(SCM(M_{p+1})=SCR(M_1)-SCR(M_{p+1})\) qui représente la réduction d’erreur quand on passe du modèle \(M_1\) au modèle \(M_{p+1}\).

TipRéécriture de la décomposition de la variance du modèle \(M_{p+1}\)

\[SCR(M_1)=(SCR(M_1)-SCR(M_{p+1}))+SCR(M_{p+1})\]

7.8 Test de \(H_0 : \mbox{ modèle } M_1\) contre \(H_1\) : modèle \(M_{p+1}\)

Statistique de test \[T_n=\frac{SCM(M_{p+1})/p}{SCR(M_{p+1})/(n-p-1)}\sim_{H_0} \mathcal{F}(p,n-p-1)\]

qui s’écrit donc aussi

\[T_n=\frac{(SCR(M_1)-SCR(M_{p+1}))/p}{SCR(M_{p+1})/(n-p-1)}\sim_{H_0} \mathcal{F}(p,n-p-1)\]

Zone de rejet

\[R_\delta = \{T_n > f_{1-\delta} \}\]

\(f_{1-\delta}\) est le quantile \(1 - \delta\) de la loi de Fisher \(\mathcal{F}(p,n-p-1)\).

Si \(t_{obs} > f_{1-\delta}\), on rejette \(H_0\) au risque \(\delta\).

7.9 Test de \(H_0 : \mbox{ modèle } M_1\) contre \(H_1\) : modèle \(M_{p+1}\)

Mise en oeuvre et conclusion :

• Si \(t_{obs} > f_{1-\delta}\), on conserve le modèle complet \(M_{p+1}\) et la contribution globale des \(p\) variables explicatives est significative : au moins une des variables explicatives contribue à expliquer la variabilité de \(Y\).
• Si \(t_{obs} < f_{1-\delta}\), on ne rejette pas \(H_0\) : le modèle \(M_1\) suffit et les variables explicatives ne contribuent pas à expliquer de manière significative \(Y\).

\(p\)-valeur :

\[p_c=\mathbb{P}_{H_0}(T_n > t_{obs})\] que l’on compare à un risque, par exemple \(\delta=5\%\)

On rejette \(H_0\) si \(p_c<\delta\).

7.10 Test avec R

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
summary(mod)

Call:
lm(formula = Weight ~ Length1 + Height + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-183.872  -17.044   -0.495   24.591  101.032 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1048.57     187.86  -5.582 4.13e-05 ***
Length1        18.96      13.26   1.430    0.172    
Height         48.50      28.31   1.713    0.106    
Width          66.72      51.73   1.290    0.215    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 63.03 on 16 degrees of freedom
Multiple R-squared:  0.808, Adjusted R-squared:  0.772 
F-statistic: 22.44 on 3 and 16 DF,  p-value: 5.624e-06

7.11 Table de l’analyse de la variance de la régression

\[M_1:\quad Y_i=\beta+ E_i, \quad E_i \overset{i.i.d.}{\sim} {\cal N}(0\,,\,\sigma^2)\]

m1 <- lm(Weight ~1, data = fish)
anova(m1, mod)

Analysis of Variance Table

Model 1: Weight ~ 1
Model 2: Weight ~ Length1 + Height + Width
  Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
1     19 331013                                  
2     16  63568  3    267445 22.439 5.624e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

7.12 Après le test

• Comme le test global \(H_0\) : \(M_1\) contre \(H_1\) : \(M_{p+1}\) est significatif, au moins une des variables explicatives contribue à expliquer \(Y\).
• Dans ce cas on peut modéliser nos données par un modèle de régression linéaire multiple (sous réserve de validation des hypothèses).
• Prendre en compte toutes les variables peut s’avérer très couteux (voir exemples en TP).

7.13 Test du modèle \(M_{q+1}\) contre le modèle \(M_{p+1}\)

WarningObjectif

• Tester si un ensemble de \(q\) variables explicatives ne suffit pas à expliquer \(Y\), avec \(q<p\).

• Autre formulation :

\[H_0\ :\ \text{modèle}\ M_{q+1}\ :Y_i=\beta+\sum_{j=1}^q\alpha_j x_{i,j}+E_i,\quad E_i\overset{iid}{\sim}{\cal N}(0,\sigma^2)\] contre l’alternative \[H_1\ :\ \text{modèle}\ M_{p+1}:Y_i=\beta+\sum_{j=1}^p\alpha_j x_{i,j}+E_i,\quad E_i\ \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

• Autre formulation :

Test de \(H_0:\ \forall j=q+1,\cdots,p,\ \alpha_j=0 \quad\) contre \(\quad H_1:\exists j=q+1,\cdots,p\ \alpha_j\neq 0\)

7.14 Description des deux modèles \(M_{q+1}\) et \(M_{p+1}\)

Modèle à \(p\) variables explicatives \(M_{p+1}\)

\[ Y_i = \beta\ +\ \sum_{j=1}^p\alpha_j\,x_{i,j}+ E_i,\quad E_i\ \overset{i.i.d.}{\sim}\ {\cal N}(0,\sigma^2).\]

Sous la forme matricielle :

\[(M_{p+1})\,: Y = X^{(p)}\theta^{(p)}\ +\ E, \quad E \sim {\cal N}(0, \sigma^2 I_n) \]

avec

\[ X^{(p)} = \begin{bmatrix} 1 & x_{1,1} & \ldots & x_{1,p}\\ 1 & x_{2,1} & \ldots & x_{2,p}\\ \vdots & \vdots & & \vdots \\ 1 & x_{n,1} & \ldots & x_{n,p} \end{bmatrix} \mbox{ et } \theta^{(p)}=[\beta,\alpha_1,\cdots,\alpha_p]^t\]

7.15 Description des deux modèles \(M_{q+1}\) et \(M_{p+1}\)

Modèle à \(q\) variables explicatives \(M_{q+1}\)

\[ Y_i = \beta\ +\ \sum_{j=1}^q\alpha_j\,x_{i,j}+ E_i,\quad E_i\ \overset{i.i.d.}{\sim}\ {\cal N}(0,\sigma^2).\]

Sous la forme matricielle :

\[(M_{p+1})\,: Y = X^{(q)}\theta^{(q)}\ +\ E, \quad E \sim {\cal N}(0, \sigma^2 I_n) \]

avec

\[ X^{(q)} = \begin{bmatrix} 1 & x_{1,1} & \ldots & x_{1,q}\\ 1 & x_{2,1} & \ldots & x_{2,q}\\ \vdots & \vdots & & \vdots \\ 1 & x_{n,1} & \ldots & x_{n,q} \end{bmatrix} \mbox{ et } \theta^{(q)}=[\beta,\alpha_1,\cdots,\alpha_q]^t\]

Le modèle à \(q\) variables explicatives \(M_{q+1}\) est emboîté dans le modèle à \(p\) variables explicatives \(M_{p+1}\) : les \(q\) variables explicatives forment un sous-ensemble des \(p\) variables explicatives.

7.16 Estimateur des paramètres dans \(M_{p+1}\) et dans \(M_{q+1}\)

Paramètres d’espérances

\[\begin{array}{ccc} \hat{\theta}^{(p)} &=& (B^{(p)},A_1^{(p)},\dots,A_{p}^{(p)})^t = \left((X^{(p)})^t\,X^{(p)}\right)^{-1}\,\!(X^{(p)})^t\,Y\\ \hat{\theta}^{(q)} & = & (B^{(q)},A_1^{(q)},\dots,A_{q}^{(q)})^t=\left((X^{(q)})^t\,X^{(q)}\right)^{-1}\,\!(X^{(q)})^t\,Y\ \end{array}\]

CautionRemarques

\(B^{(p)}\neq B^{(q)}, A_1^{(p)}\neq A_1^{(q)},\cdots\)

7.17 Prévision et résidu dans \(M_{p+1}\) et dans \(M_{q+1}\)

NotePrévision dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[\begin{align}\widehat{Y_i}(M_{p+1})&=B^{(p)}+\sum_{j=1}^p A^{(p)}_jx_{i,j}\\ \widehat{Y_i}(M_{q+1})&=B^{(q)}+\sum_{j=1}^q A^{(q)}_jx_{i,j} \end{align}\]

NoteRésidu dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[\begin{align}E_i(M_{p+1})&=Y_i-\widehat{Y_i}(M_{p+1})=Y_i-(B^{(p)}+\sum_{j=1}^p A^{(p)}_jx_{i,j})\\ E_i(M_{q+1})&=Y_i-\widehat{Y_i}(M_{q+1})=Y_i-(B^{(q)}+\sum_{j=1}^q A^{(q)}_jx_{i,j}) \end{align}\]

7.18 SCR dans \(M_{p+1}\) et dans \(M_{q+1}\)

NoteSCR dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[SCR(M_{p+1}) = \sum_{i=1}^{n}\left(Y_i-B^{(p)}-\sum_{j=1}^{p}A^{(p)}_j\,x_{i,j}\right)^2\ =\ ||Y - X^{(p)}\hat{\theta}^{(p)}||^2\]

\[SCR(M_{q+1}) = \sum_{i=1}^{n}\left(Y_i-B^{(q)}-\sum_{j=1}^{q}A^{(q)}_j\,x_{i,j}\right)^2\ =\ ||Y - X^{(q)}\hat{\theta}^{(q)}||^2\]

NoteEstimation de \(\sigma^2\) dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[S^2_{(M_{p+1})}=\frac{SCR(M_{p+1})}{n-p-1}\\ S^2_{(M_{q+1})}=\frac{SCR(M_{q+1})}{n-q-1}\]

7.19 Résultat important

TipThéorème important du cours

Sous les hypothèses du modèle linéaire Gaussien multiple, et si on définit \(SCE=SCR(M_{q+1})-SCR(M_{p+1})\) alors

• \(\frac{SCR(M_{p+1})}{\sigma^2}\, \sim\,\chi^2_{n-p-1}\)
• Sous le modèle \(M_{q+1}\), c’est à dire sous \(H_0\) :

\[\frac{SCR(M_{q+1})}{\sigma^2} \underset{H_0}{\sim}\,\chi^2_{n-q-1}\quad;\quad \frac{SCE}{\sigma^2} \underset{H_0}{\sim}\,\chi^2_{p-q}\\\]

et \(SCE\) et \(SCR(M_{p+1})\) sont indépendants ce qui permet de déduire que

\[\frac{\left(SCR(M_{q+1}) - SCR(M_{p+1})\right)\Bigl/(p-q)}{SCR(M_{p+1})\Bigl/~(n-p-1)} \sim_{H_0}^{}~\mathcal{F}(p-q\,,\,n-p-1)\]

7.20 Test de \(H_0\) : modèle \(M_{q+1}\) contre \(H_1\) : modèle \(M_{p+1}\)

TipStatistique de test

\[ {T_n=\frac{(SCR(M_{q+1})-SCR(M_{p+1}))/(p-q)}{SCR(M_{p+1})/(n-p-1)}\sim_{H_0} \mathcal{F}(p-q,n-p-1)} \] - \(SCR(M_{q+1})-SCR(M_{p+1})\) représente la réduction d’erreur quand on passe de \(M_{q+1}\) à \(M_{p+1}\).

TipZone de rejet

\[R_\delta = \{T_n > f_{1-\delta} \} \] où \(f_\delta\) est le quantile d’ordre \(1-\delta\) de la \(\mathcal{F}(p-q,n-p-1)\).

On rejette \(H_0\) si \(t_{n}> f_{1-\delta}\)

7.21 Test de \(H_0\) : modèle \(M_{q+1}\) contre \(H_1\) : modèle \(M_{p+1}\)

TipInterprétation

Si \(t_n > f_{1-\delta}\), on conserve le modèle complet \(M_{p+1}\), on considère que le passage de \(M_{q+1}\) à \(M_{p+1}\) est significatif : au moins une variable explicative parmi \(x_{q+1},\cdots,x_p\) a une influence significative sur \(Y\) (en plus de \(x_1,\cdots,x_q\)).

Si \(t_n \leq f_{1-\delta}\), on ne rejette pas \(H_0\), on considère que le passage de \(M_{q+1}\) à \(M_{p+1}\) n’est pas significatif : l’influence des variables explicatives \(x_{q+1},\cdots,x_p\) n’est pas significative pour expliquer \(Y\).

TipCalcul de la \(p\)-valeur

\[p_c=\mathbb{P}_{H_0}(T_n > t_n)=\mathbb{P}(\mathcal{F}(p-q,n-p-1)>t_n)\]

7.22 Table de l’analyse de la variance de la régression

Avec R, commande anova

\[\begin{array}{|c|c|c|c|c|c|c|} \hline Source & ddl & SCR& ddl & SCE\quad expliquée& Statistique & p-value\\ & résiduelle& && par\quad M_{p+1} & de\quad test & \\ \hline \text{Modèle} & n-q-1 & SCR(M_{q+1}) & &&&\\ M_{q+1} &&&&&&\\ \hline \text{Modèle} & n-p-1 & SCR(M_{p+1}) & p-q & SCE=SCR(M_{q+1}) & t_n=&\mathbb{P}(\mathcal{F}(p-q,n-p-1)>t_n)\\ M_{p+1} &&&& - SCR(M_{p+1}) &\frac{SCE/(p-q)}{SCR(M_{p+1})/(n-p-1)} &\\ \hline \end{array}\]

7.23 Retour à l’exemple

On va tester

\[\begin{align} &M_4 : Y_i=\beta+ \alpha_1x_{i,1}+\alpha_2x_{i,2}+\alpha_3x_{i,3}+ E_i \\ &M_3 : Y_i=\beta+ \alpha_1x_{i,1}+\alpha_2x_{i,2} +E_i \end{align}\]

On va définir le modèle \(M_3\)

m3 <- lm(Weight ~ Length1 + Width, data = fish)
summary(m3)

Call:
lm(formula = Weight ~ Length1 + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-211.090  -14.742    2.901   25.347  104.568 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -998.742    195.859  -5.099 8.91e-05 ***
Length1       35.713      9.449   3.779   0.0015 ** 
Width         97.835     51.111   1.914   0.0726 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 66.52 on 17 degrees of freedom
Multiple R-squared:  0.7728,    Adjusted R-squared:  0.746 
F-statistic:  28.9 on 2 and 17 DF,  p-value: 3.391e-06

7.24 Sorties R

anova(m3,  mod)

Analysis of Variance Table

Model 1: Weight ~ Length1 + Width
Model 2: Weight ~ Length1 + Height + Width
  Res.Df   RSS Df Sum of Sq      F Pr(>F)
1     17 75223                           
2     16 63568  1     11655 2.9335 0.1061

anova(m1,  m3)

Analysis of Variance Table

Model 1: Weight ~ 1
Model 2: Weight ~ Length1 + Width
  Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
1     19 331013                                  
2     17  75223  2    255790 28.904 3.391e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

7.25 Remarques sur les modèles emboités

CautionRemarques

• On peut toujours comparer des modèles emboités l’un dans l’autre.
• On peut donc comparer tous les modèles au modèle complet (structure la plus riche).
• On peut de même comparer tous les modèles au modèle constant (structure la moins riche).
• On ne peut pas comparer des modèles non emboîtés à l’aide de ce test :